# This code replicates Table 10
rm(list=ls())
setwd("/Users/anselmrink/Dropbox/Eigener Research/Peru/Authority/Analysis/Quantitative data/00 Repl Arch")
library(foreign)
data <- read.dta("data.dta")

# Collapse data to village level
data <- data[c("village", "treatment_catholic", "treatment_adventist", "treatment_mixed", "treatment_peruana", "treatment_maranatha", "total_pop",	"male_census",	"female_census",	"x1_below",	"x1_4",	"x5_14",	"x15_64",	"x65plus",	"natives_census",	"migrants_census",	"foreigners",	"handycap",	"blind",	"mental_issues",	"polio",	"underdeveloped_extremities",	"overdeveloped_extremities",	"other_handycaps",	"analphabets_total",	"analphabets_male",	"analphabets_female",	"no_educ",	"preschool",	"primary_educ",	"secondary_educ",	"higher_educ",	"workers6_14",	"workers_15plus",	"workers_15plus_employed",	"workers_15plus_unemployed",	"occ_agriculture",	"occ_manu_mine_constr",	"occ_retail",	"occ_hawker",	"occ_other_services",	"occ_other",	"cat_occ_wage_earner",	"cat_occ_ind",	"cat_occ_employer",	"cat_occ_fam",	"cat_occ_household",	"econ_primary",	"econ_secondary",	"econ_tertiary",	"living_together",	"married_census",	"alone",	"other",	"head_household",	"head_household_male",	"head_household_female",	"average_children_per_mother",	"women_with_more_4_children",	"single_mothers",	"single_mothers_12_19",	"single_mothers_20_29",	"single_mothers_30_49",	"young_mothers",	"flats_with_unoccupied_homes",	"households",	"average_household_size",	"ind_house",	"makeshift",	"other_1",	"ownership_own",	"ownership_rental",	"ownership_occupied",	"ownership_other",	"walls_cement",	"walls_quincho",	"walls_stone",	"walls_wood",	"walls_mat",	"walls_other",	"roof_concrete",	"roof_corrugated_steel",	"roof_mats",	"roof_straw",	"roof_other",	"water_private",	"water_public_well",	"water_tankwagon",	"water_other",	"sanitary_connected_grid",	"sanitary_privy_pit",	"sanitory_other",	"sanitory_no",	"flats_with_electricity",	"flats_no_electricity",	"flats_with_1_room",	"households_only_dormitories",	"households_share_toilet",	"households_with_commerical_roms",	"no_domestic_appliances",	"only_radio",	"radio_tv",	"sewing_machine",	"fridge",	"tricycle",	"with_4plus_appliances")]
data <- aggregate(data, by=list(data$village), FUN=mean, na.rm=TRUE)

## Define variables
treats <- c("treatment_catholic", "treatment_adventist", "treatment_mixed", "treatment_peruana", "treatment_maranatha")
vars <- c("total_pop",	"male_census",	"female_census",	"x1_below",	"x1_4",	"x5_14",	"x15_64",	"x65plus",	"natives_census",	"migrants_census",	"foreigners",	"handycap",	"blind",	"mental_issues",	"polio",	"underdeveloped_extremities",	"overdeveloped_extremities",	"other_handycaps",	"analphabets_total",	"analphabets_male",	"analphabets_female",	"no_educ",	"preschool",	"primary_educ",	"secondary_educ",	"higher_educ",	"workers6_14",	"workers_15plus",	"workers_15plus_employed",	"workers_15plus_unemployed",	"occ_agriculture",	"occ_manu_mine_constr",	"occ_retail",	"occ_hawker",	"occ_other_services",	"occ_other",	"cat_occ_wage_earner",	"cat_occ_ind",	"cat_occ_employer",	"cat_occ_fam",	"cat_occ_household",	"econ_primary",	"econ_secondary",	"econ_tertiary",	"living_together",	"married_census",	"alone",	"other",	"head_household",	"head_household_male",	"head_household_female",	"average_children_per_mother",	"women_with_more_4_children",	"single_mothers",	"single_mothers_12_19",	"single_mothers_20_29",	"single_mothers_30_49",	"young_mothers",	"flats_with_unoccupied_homes",	"households",	"average_household_size",	"ind_house",	"makeshift",	"other_1",	"ownership_own",	"ownership_rental",	"ownership_occupied",	"ownership_other",	"walls_cement",	"walls_quincho",	"walls_stone",	"walls_wood",	"walls_mat",	"walls_other",	"roof_concrete",	"roof_corrugated_steel",	"roof_mats",	"roof_straw",	"roof_other",	"water_private",	"water_public_well",	"water_tankwagon",	"water_other",	"sanitary_connected_grid",	"sanitary_privy_pit",	"sanitory_other",	"sanitory_no",	"flats_with_electricity",	"flats_no_electricity",	"flats_with_1_room",	"households_only_dormitories",	"households_share_toilet",	"households_with_commerical_roms",	"no_domestic_appliances",	"only_radio",	"radio_tv",	"sewing_machine",	"fridge",	"tricycle",	"with_4plus_appliances")

## Write helper function
get_stats <- function(var, treat, data) {
  c(mean(data[, var][data[, treat] ==1 ], na.rm = T), sd(data[, var][data[, treat] ==1 ], na.rm = T))
}

## Loop through variables
btable <- lapply(treats, function(y) t(sapply(vars, function(x) {
  get_stats(x, y, data = data)
})))

## To data frame
btable <- data.frame(do.call(cbind, btable))

## Correct percentages
perc <- c("male_census",	"female_census",	"x1_below",	"x1_4",	"x5_14",	"x15_64",	"x65plus",	"natives_census",	"migrants_census",	"foreigners",	"handycap",	"blind",	"mental_issues",	"polio",	"underdeveloped_extremities",	"overdeveloped_extremities",	"other_handycaps",	"analphabets_total",	"analphabets_male",	"analphabets_female",	"no_educ",	"preschool",	"primary_educ",	"secondary_educ",	"higher_educ",	"workers6_14",	"workers_15plus",	"workers_15plus_employed",	"workers_15plus_unemployed",	"occ_agriculture",	"occ_manu_mine_constr",	"occ_retail",	"occ_hawker",	"occ_other_services",	"occ_other",	"cat_occ_wage_earner",	"cat_occ_ind",	"cat_occ_employer",	"cat_occ_fam",	"cat_occ_household",	"econ_primary",	"econ_secondary",	"econ_tertiary",	"living_together",	"married_census",	"alone",	"other"                     ,	 "head_household_male",	"head_household_female",	"women_with_more_4_children",	"single_mothers",	"single_mothers_12_19",	"single_mothers_20_29",	"single_mothers_30_49",	"young_mothers",	"makeshift",	"other_1",	"ownership_own",	"ownership_rental",	"ownership_occupied",	"ownership_other",	"walls_cement",	"walls_quincho",	"walls_stone",	"walls_wood",	"walls_mat",	"walls_other",	"roof_concrete",	"roof_corrugated_steel",	"roof_mats",	"roof_straw",	"roof_other",	"water_private",	"water_public_well",	"water_tankwagon",	"water_other",	"sanitary_connected_grid",	"sanitary_privy_pit",	"sanitory_other",	"sanitory_no",	"flats_with_electricity",	"flats_no_electricity",	"flats_with_1_room",	"households_only_dormitories",	"households_share_toilet",	"households_with_commerical_roms",	"no_domestic_appliances",	"only_radio",	"radio_tv",	"sewing_machine",	"fridge",	"tricycle",	"with_4plus_appliances")
btable[rownames(btable) %in% perc, ] <- btable[rownames(btable) %in% perc, ] * 100

## Add varnames
btable$variable <- rownames(btable)

## Add N
n <- sapply(vars, function(x) sum(!is.na(data[, x])))
btable <- cbind(n, btable)

## Change colnames
rownames(btable) <- c("Total population",	"Male",	"Female",	"Age below 1",	"Age 1-4",	"Age 5-14",	"Age 15-64",	"Age over 64",	"Natives",	"Migrants",	"Foreigners",	"Handycap",	"Blind",	"Mental issues",	"Polio",	"Underdeveloped extremities",	"Overdeveloped extremities",	"Other handycaps",	"Analphabets (total)",	"Analphabets (male)",	"Analphabets (female)",	"No education",	"Preschool",	"Primary education",	"Secondary education",	"Higher education",	"Workers 6-14",	"Workers over 14",	"Workers over 14 (employed)",	"Workers over 14 (unemployed)",	"Occupation: Agriculture",	"Occupation: Construction",	"Occupation: Retail",	"Occupation: Hawker",	"Occupation: Other services",	"Occupation: Other",	"Wage earner",	"Independent worker",	"Employer",	"Family Employment",	"Working in household",	"Employed in primary sector",	"Employed in secondary sector",	"Employed in tertiary sector",	"Living together",	"Married",	"Living alone",	"Living in other form",	"Heads of household",	"Head household (male)",	"Head household (female)",	"Average children per mother",	"Women over 4 children",	"Single mothers",	"Single mothers 12 19",	"Single mothers 20-29",	"Single mothers 30-49",	"Young mothers",	"Flats incl. unoccupied homes",	"Households",	"Average household size",	"Independent houses",	"Makeshift",	"Other house",	"Ownership: own",	"Ownership: rental",	"Ownership: occupied",	"Ownership: other",	"Walls: cement",	"Walls: quincho",	"Walls: stone",	"Walls: wood",	"Walls: mat",	"Walls: other",	"Roof: concrete",	"Roof: corrugated steel",	"Roof: mats",	"Roof: straw",	"Roof: other",	"Water: private",	"Water: public well",	"Water: tankwagon",	"Water: other",	"Sanitary: connected grid",	"Sanitary: privy pit",	"Sanitory: other",	"Sanitory: none",	"Flats with electricity",	"Flats no electricity",	"Flats with only 1 room",	"Households only dormitories",	"Households with shared toilet",	"Households with comm rooms",	"No domestic appliances",	"Only radio",	"Radio and TV",	"Sewing machine",	"Fridge",	"Tricycle",	"More than 4 appliances")
cn <- paste0(rep(treats, each = 2), rep(c("_mean", "_sd"), 2))
colnames(btable)[c(-1, -ncol(btable))] <- cn
btable <- data.frame(round(btable[, -which(colnames(btable) == "variable")], 1))
btable

library(xtable)
xtable(btable, digits=1)
